Exact results for spatial decay of correlations in low-dimensional insulators II 
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We study decay rates of one-body reduced density matrices in insulators, described by a tight- 
binding model, where not only an external potential but also hoppings are spatially modulated. We 
determine analytically the power in the power law and the correlation length in D = 1 case and 
in several lattice directions in D — 2 case. Unlike the uniform hopping case, in D = 1 system and 
in some directions of D — 2 system the correlation length is not determined uniquely by the gap. 
Moreover, a crossover from D — 2-decay rates to D — 1 ones is investigated. 
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The localization properties of one-body reduced den- 
sity matrices (DM) are of vital importance for rapidly de- 
veloping computational techniques, used for calculating 
various characteristics of crystalline solidsiS^. But the 
interest in large-distance decay rates of DM dates from 
1959, when Kohn^ determined an exponential factor in 
the large-distance asymptotics of DM in one-dimensional 
(D = 1) crystals. Recent works concerning this problem 
include R e fsiii£iL2i2ii2iiiii£. In particular, another major 
step, after Kohn's result, was made by He and Vander- 
bilt who discovered a power-law factor in the dacay rate 
of D = 1 crystals. Extensions to systems in dimensions 
D > 1 have been obtained by Des Cloizeausi^, Taraskin 
et aliSiii, and in our recent paper—, where we have shown 
to what extent the large-distance asymptotics of DM is 
sensitive to lattice directions in D = 2 crystals. 

In all the papers quoted above, the pattern of hoppings 
is spatially uniform. Here, we consider a system of spin- 
less fermions on a lattice whose hoppings are spatially 
modulated, subjected to a modulated external potential. 
The Hamiltonian of the considered system reads: 



H = ) Uiafai + ^2 Uj (afaj + h.c. 
* (i,3) 



(1) 



In the above expression, i,j represent the lattice sites of 
a D-dimensional simple cubic lattice, stand for the 

pairs of nearest neighbors on this lattice, and in the sum 
over nearest neighbors each pair is counted once. The 
operators af (a^) create (annihilate) a spinless fermion in 
a single-particle orbital \i), belonging to an orthonormal 
basis. As j runs through the underlying lattice, the value 
of the external potential at site j is modulated with the 
amplitude u, 



Uj = U+(-l)^u, 



(2) 



where for a site j = (ji,...,j D ), Oj = J2i=i \jl\ is the 
(noneuclidean) distance of this site from the origin. Since 
a D-dimensional simple cubic lattice consists of two in- 
terpenetrating sublattices: the even one - with aj even, 
and the odd one - with aj odd, the external potential 
is uniform at each sublattice. Its value changes by ±2it 
on passing from one sublattice to the other one, and the 



average over the whole lattice is U . Similarly, the sym- 
metric hopping matrix tij is modulated: 



t ii±ei =t[l±(-l) CTi r] 



(3) 



with the mean t and the amplitude tr, where ej, I = 
1, D, is the unit vector with only one nonzero compo- 
nent - the l-th. one. Without any loss of generality, we 
can suppose that u > 0, t > 0, and 1 > r > 0. 

The considered model can describe electronic prop- 
erties of crystal surfaces, possibly containing adsorbed 
atoms. The pattern of hoppings © may arise when 
atoms form a distorted square lattice, that is one of the 
sublattices is shifted with respect to the other. 

Under the periodic boundary conditions, the Hamilto- 
nian is block-diagonalizcd by the plane wave orbitals, 
|fe), with the wave vector in the first Brillouin zone of the 
lattice. Specifically, shifting zero of the energy scale to 
U and expressing all the energies involved in the units of 
the average transfer integral t, we obtain the upper, Aj£, 
and the lower, A^, bands of cigcnencrgics: 



A 



± 



(4) 



In Eq.{3J, the wave vector k is restricted to the first 
Brillouin zone of one of the sublattices (so that there 
is as many eigenvalues as sites in the lattice), — 

2J2\Zi cosfc;, (k = 2r J2\Zi sinfc;. The two bands arc 
mutually symmetric about zero and are separated by a 
gap go- If the system given by QJ is half-filled and 
gD > 0, then it is an insulator in the sense of the standard 
band theory. The normalized eigenvectors correspond- 
ing to the eigenvalues Aj7 are linear combinations of the 
vectors |fc) and \k + n), where tt stands for the vector 
with all the components equal to it. By means of these 
eigenvectors, one can calculate the zero-temperature non- 
diagonal elements of DM for the half-filled system, which 
in the thermodynamic limit read 

1 D 



(afa i+r )=--(-iy i n D (r), if a. 



if av. = 2m+ 1, and 

2m, (5) 



2 



where 



n D (r) =2(2%)-° [ dkexp (ikr) A k \ 

J B.Z. 



(6) 



and S{R(r) = TZ(n + l,...,r D )+ H(r x - 1, ... , r D ), 
V{R{r) = lZ(ri + l, . . . , r\o) — TZ(n — 1, . . . , ro), and so on 
for I = 2, . . . , D. In Eq. ©, the D-dimensional integral 
is taken over the first Brillouin zone of a sublattice and 
r /: 0. At odd sites (oy = 2m + 1), the function TZd 
vanishes. 

At the level of the Hamiltonian it is clear that 
for r = we are back in the the uniform hopping case 
of Refii£. Note however, that here is twice as large 
as in Ref4£, and consequently an extra factor 2 appears 
in the definition of 7?._D(r). On the other hand, if r = 
1, the system looses its connectivity with respect to the 
hopping. Specifically, for D = 1, we obtain a system of 
independent dimers whose correlation functions vanish at 
distances greater than one, while for D = 2 - a system 
of independent quasi-one-dimensional zig-zag chains. 

In what follows, we shall investigate the large-distance 
decay rate of the function IZoir) in D = 1,2 systems. 

It is instructive to analyze first the one-dimensional 
case. The gap in the spectrum is gi = ASi, with 
Sf = p\ + r 2 , and Pi = u/2. Apparently, the role of 
the rescaled modulation amplitudes, u and r, in creating 
the gap is symmetric. Moreover, g\ > unless u = r = 0. 
For later purposes, it is convenient to introduce, besides 
the rescaled gap Si, one more parameter describing the 
spectrum: for the lower band this parameter amounts to 
the ratio, Ki, of the upper edge of the band to the lower 



edge, i.e. 



(1 - t 2 )/(1 +(3{). The pair (<5i,ki) can 



be thought of as new coordinates, equivalent to the co- 
ordinates (Pi,t). The specific combination of Pi and t, 
denoted Ki, appears naturally in an expression for the 
correlation function (JSJ. Provided that Ki ^ 0, the func- 
tion lZi(r) at even points, r = 2m, can be written as 



Ki(2m) = 



1 



1 



(-1)™ 1 



dk 



lm-l/2 



cos(2?nfc) 
1 — k\ sin 2 k 



where T stands for the Euler T-function, F - for the 
Gauss hypergeometric function, and Q v - for the Legen- 
dre function of second kincU^. Eq. implies the follow- 
ing large-m asymptotics: 



■Ri(2m) 



(-l) m (1 - k 2 ) 1 / 4 exp(-2m/£) 



./2tt Si ^2m 
where the inverse correlation length £ -1 is given by 

1 



r 



ln 



Kl 



(8) 



(9) 



We note that the correlation function (afai+r) is con- 
tinuous at r = t = 1, and r = j3\ = 0. Out of the 



three above cases, the asymptotics (0 is continuous only 
at t = 0, where we recover the asymptotics obtained in 
Re02. 

The new variables reveal clearly the nature 

of the relation between £, Ki and 5\. It follows from 
© that Ki determines uniquely £ and vice versa. But 
the relation between £ (or ki) and <5i is a multi- valued 
one. Specifically, to given Si there corresponds an inter- 
v al [Kf'" (Ji),K7 aa (Ji)] of ki values, where K^ m {Si) = 
s/l - Si if Si < 1, and K\ nm (Si) = if Si > 1, and 
K™ a:r (<5i) = + S\. Consequently, the inverse corre- 

lation length £ _1 can assume any value from the interval 
[f-!(Kf f-^ftf *"(£!))]. As Si tends to zero, the 

edges of the interval tend to zero as well but its width 
shrinks even faster, what results in a fairly good local- 
ization of £ , for small Si- Specifically, expanding the 
bounds at Si = 0, 



C 1 (KT ax (Si))^Si--Sl + -8l + 



40 



(10) 



(11) 



we sec that to the first order in Si the bounds coincide and 
the indeterminacy-interval width is Sf/2, in the leading 
order. 

For Si > 1, £ -1 can assume any value not less than 
£ _1 (K5 na:E (i5i)), and this lower bound blows up as Si — > 
oo. 

The described above indeterminacy of vanishes if 
one of the parameters, Pi or r, is kept constant, since 
such conditions turn the relations between «i and into 
single valued functions. Of particular interest is the be- 
havior of in two such situations, motivated physi- 
cally. Firstly, a weak modulated external potential (/3i 
varies in a vicinity of zero) is imposed on a connected 
system with modulated hopping (0 < r < 1, fixed). Sec- 
ondly, a weak modulated hopping is imposed (r varies 
in a vicinity of zero) on a system in a modulated poten- 
tial (0 < Pi, fixed). In the both cases, Si is separated 
1 ) , (7) from zero. Therefore, to describe the behavior of it 
is more convenient to consider the expansions in Pi or r, 
respectively. Specifically, in the first situation 



In 



1 



IPI 
t 2 



1 



Hi 



while in the second one 



ln(A + ^l + Pl) + 



y/TT^T 2 



(12) 



(13) 



Clearly, in the above expansions it is not possible to 
go to zero with the fixed parameters. To obtain the 
corresponding expansions about zero, we note that the 
function n™ ax (5i) = Ki(Pi,t)\/3 1= Si,t=o- Therefore, for 
r = the expansion in Pi coincides with the expansion 
(fTUft of the lower bound for Similarly, k™"(<5i) = 



3 



T)\p 1 —o tT —s 1 , hence for /3i = the expansion in r 
coincides with the expansion of the upper bound for 

r 1 . 

Yet another interesting regime is, when a connected 
system (t 7^ 1) is subjected to a strongly modulated ex- 
ternal potential. Then, the large j3\ asymptotics of £ _1 
reads: 



r 



In 



In ft 



1 



4/3? 



(14) 



The two-dimensional system differs from the one- 
dimensional one in many respects, due to different ge- 
ometry and connectivity of the underlying lattice. When 
hopping is modulated, an important difference appears 
already in the spectrum: the gap g 2 = 8(3 2 , /?2 = u/A, is 
independent of the amplitude r, and vanishes unless the 
external potential is modulated. 

The function lZ 2 (r) can be cast into the form 



1 f 2 f 2 dxdy cos(2mcc) cos(2ny) 
Wo Jo y/(3l + [t 2 + ( 1 - T 2 ) cos 2 y] cos 2 



= , (15) 



where we introduced the integers m and n such that 
2m = r\ — r 2 , 2n = r\ + r 2 . For r = 0, 712(1") co- 
incides with that of Refw On the other hand, in the 
case of quasi-one-dimensional zig-zag chains (i.e. r = 1), 
lZ 2 (r) is nonzero only along direction n — 0, and in this 
direction it amounts to 72i(r) but with twice as large 
hopping as in the uniform one-dimensional case. 

By expanding the square root in the denominator, the 
double integral of (|15|l can be related to a single integral 
involving a hypergeometric function ^F^^: 



^2 (r) 



2\-n 



T(m + 1/2) 



X3F2 |m + -,TO- 
for to > n, and 



2tt (2/3 2 ) 2m+1 r(n+l/2)r(m-n+l) 

Jt 2 

1 1 



-, m + 1; 2m+l, to— n + 1; 



v 



(-1)"(1-T 2 )- 



r(2n + l) 



2tt (2/? 2 ) 2 " +1 r(n+m+l)r(n-m+l) 



d«(l-«) n - 1 / 2 («-T 2 ) n - 1 / 2 



n+^,n+^,n+l;n + m+l,n-m+l;--^]. 
for n > m. 

Using representations (|16[) . I|17|l . we have derived large- 
distance asymptotic behavior of 72,2 (r) in diagonal direc- 
tions of a connected system (r 7^ 1). Specifically, in 
direction n = it follows from (jlf)|) that 



72.2(m, — to) 



(-1)' 



cxp(— 2m/£) 



2tt/3 2 yf— ^2 



2m 



(18) 




FIG. 1: The inverse correlation length of D=2 system in direc- 
tion n = m against u, for different values of r. Symbols rep- 
resent the curves obtained from fitting numerically lZ(2m, 0). 
The continuous lines are derived from approximate formula 
0, the upper one for r = 0.25, the lower one for r = 0.05. 
The dashed curve is derived from an exact formula for r = 0. 
In the inset we display the corresponding log-log plots. 



C L = M\i 



■Pi)- 



(19) 



Apparently, in direction n — the gap determines £. The 
above asymptotics differs only by the factor (1 — t 2 ) -1 / 2 
from that in the uniform hopping case^, hence it is con- 
tinuous at r = 0. More interesting is a neighborhood of 
t = 1, where we would like to observe how our D = 2 
system approaches the corresponding D = 1 system. In 
particular, an interesting question is how the power law 
r 7 , with r = 2m and 7 = 1 in l|18|) . transforms into \j2m 
in as r — ► 1? Our numerical calculations suggest 
that 7 behaves discontinuously: it grows to some value 
above 1, and then jumps to 1/2 at r = 1. 
Then, in direction m = formula 117|) gives 



722(n, n) 



(-1)" exp(-2n/Q 
27T/? 2 2n 



(20) 



where £ _1 is given by 10 with n\ replaced by k 2 and k 2 
obtains from k\ by changing f3\ into (3 2 . Similarly, we 
introduce S 2 , <5| = /?| + r 2 . While 5 2 and K2 lack the 
physical interpretation of Si and K\, the formal relation 
between them is fairly analogous to that between the 
latter parameters. Consequently, in direction to = the 
gap does not determine 

In directions other than diagonal ones, we have not 
been able to establish the asymptotic behavior of the 
correlation function. But with the help of numerical cal- 
culations, we have obtained an approximate asymptotics 



4 



in the axial direction m = n. This asymptotics can be 
considered as an upper bound for the true one, and it 
approaches the true one as u grows and/or r decreases. 
Specifically, after expanding the hypergeometric function 
of H6|) (or (|17[l ) in a series 1 ^, and carrying out the in- 
tegration, we get a power series in — 1//3| whose coeffi- 
cient by power j depends on the hypergeometric function 
F(—j, m+1/2; 2m+l; 1 — t 2 ). One can verify numerically 
that if t <C 1 and j <C m, then this function is approx- 
imately equal to F (-j, m + 1/2; 2m + l; 1) (l + t 2 ) j . On 
substituting this approximation to the series and sum- 
ming up the series, we arrive at an expression similar to 
formula (9) of Refpi^, which in turn constitutes an ap- 
proximation to 7?.2(2m, 0). The obtained approximation 
for the correlation function is the better the smaller is 
r and/or the larger is u (see Fig^J. In particular, if 
t = 0, it becomes exact for any fa- On the basis of the 
results of Ref^, the approximate asymptotic behavior of 
72-2 (2m., 0) has the form: 



ft 2 (2m,0)! 



(-1)' 



1 



cxp(— 2m/£) 



irVWi $1 + ^ 2m 

a — 1 7r 
a+I~4 



x cos 771 arccos ■ 



(21) 



where 



C 1 = In 



(y/a+Vl + a) 



■ft) 



(22) 



As can be seen in Fig^ the above approximate £ _1 is a 
lower bound for the numerically exact one. The asymp- 
totic behaviors of the inverse correlation length ilL'L'l) , for 



small and large gaps, read: 



dn 



1- 



ln- 



1-T 2 

4 



1 P2 



vT 



In (3 2 + 



16/3, 



35 r 4 
512$ 



,(23) 



respectively. 

To summarize, continuing our studies of the decay 
rates of one-body reduced density matrices in insulators, 
by means of analytic and numerical methods, we have 
investigated a tight-binding model, where not only an ex- 
ternal potential but also hoppings can be spatially modu- 
lated. The general result of our investigations, presented 
in this report and in Ref^, reads: the form of the large- 
distance decay is always a power law times an exponen- 
tial, unless the hoppings and the external potential are si- 
multaneously uniform. We have determined analytically 
the power in the power law and the correlation length for 
D = 1 system and in several directions for D = 2 sys- 
tem. In comparison to the uniform hopping case, the con- 
sidered here system with modulated hoppings exhibits 
interesting novel features. Firstly, since the correlation 
length depends on two energy parameters (controlling 
the amplitudes of spatial modulation of the external po- 
tential and the hopping), in D = 1 case and in some 
directions of D = 2 case the gap does not determine the 
correlation length. Secondly, by adjusting the hopping- 
modulation amplitude suitably, the D = 2-system can 
be decoupled into quasi-onc-dimcnsional zig-zag chains, 
which has made possible investigating a crossover from 
D = 2-decay rates to D = 1 ones. 
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